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Inverse bicontinuous cubic phases are ubiquitous in lipid-water mixtures and consist of a lipid 
bilayer forming a cubic minimal surface, thereby dividing space into two cubic networks of water 
' channels. For small hydrocarbon chain lengths, the monolayers can be modeled as parallel surfaces 

, to a minimal midsurface. The bending energy of the cubic phases is determined by the distribution 

of Gaussian curvature over the minimal midsurfaces which we calculate for seven different structures 
(G, D, P, I-WP, C(P), S and F-RD). We show that the free-energy densities of the structures G, 
D and P are considerably lower than those of the other investigated structures due to their narrow 
distribution of Gaussian curvature. The Bonnet transformation between G, D, and P implies that 
these phases coexist along a triple line, which also includes an excess water phase. Our model 
includes thermal membrane undulations. Our qualitative predictions remain unchanged when higher 
order terms in the curvature energy are included. Calculated phase diagrams agree well with the 
experimental results for 2:1 lauric acid/dilauroyl phosphatidylcholine and water. 
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I. INTRODUCTION 



When dissolved in aqueous solvent , lipids self- assemble into an amazing variety of different structures as a function 
I of concentration and temperature . Most prominent is the lamellar phase, which consists of a stack of lipid bilayers 
(~| . separated by layers of water. In fact, the lipid bilayer constitutes the basic building block of biological membranes; in 
O ' eucaryotic cells, it not only envelopes the cell itself, but also separates its nucleus and its organelles from the cytosol, 
, and builds up the endoplasmic reticulum, the Golgi apparatus and the uni-lamellar vesicles which serve as reaction 
chambers and transport vehicles. The predominance of the lamellar phase at ambient temperatures derives from the 
I I fact that in contrast to surfactants, which usually have large headgroups and form micelles, lipids have rather bulky 
^ ' hydrocarbon chains. A so-called bilayer lipid has a chain region cross section which is of similar size as the area per 
^ ' head group; thus the spontaneous curvature is small and the lipid forms a lamellar phase. 

I If temperature is increased, the size of the head groups' hydration shell decreases, spontaneous curvature increases 
' and the lipid monolayers tend to curve towards the water regions. Since the tendency of lipid monolayers to curve is 
, frustrated in the lamellar phase, it is often found to transform into a inverse hexagonal phase which consists of a two- 
■ dimensional array of inverse cylindrical micelles. Spontaneous curvature can also be increased by changing molecular 
, architecture, by adding lipids with bulkier chains or by increasing salt concentration (since salt screens electrostatic 
"Y^ repulsion between charged headgroups) [Q. Lipids with spontaneous curvature are often called non-bilayer lipids. For 
^ large spontaneous curvature, a cubic phase of inverse spherical micelles is sometimes observed The inverse phases 
have been denoted type II (water-in-oil) phases in contrast to type I (oil-in-water) phases which dominate surfactant 
' phase behavior ||^. Although there are some lipids which also form type I phases, we focus here on type II phases 
since these are the ones formed by the majority of biologically relevant lipids. 

One of the intriguing aspects of the polymorphism of lipid-water mixtures is the existence of yet another structural 
^ ' type, which often occurs between the lamellar and the inverse hexagonal phase. Like the lamellar phase, these inverse 
bicontinuous cubic phases (IBCPs) consist of lipid bilayers. However, now a single lipid bilayer extends throughout 
the whole of space, and divides it into two disconnected but interpenetrating labyrinths of cubic symmetry, which 
both are filled with water. Thus the structure can also be considered to be built from water channels (that is inverse 
?H cylindrical micelles) as is the inverse hexagonal phase, only that now the channels meet in vertices with a coordination 
number related to the topology of the cubic phase. Although the IBCPs unite the structural elements of the lamellar 
and the inverse hexagonal phases, in contrast to them they are optically isotropic and highly viscous. Like the inverse 
hexagonal phase, the IBCPs are stabilized by spontaneous curvature. Therefore they are formed by non-bilayer lipids 
and their stability can be controlled by changing the lipids' spontaneous curvature. 

The generic occurrence of IBCPs in lipid-water mixtures and their property of dividing space into two interwoven 
aqueous compartments lends itself readily to speculations on their biological relevance, e.g. on their possible occurrence 
in the endoplasmic reticulum or the Golgi apparatus |^j^LA recent survey of a vast number of TEM pictures suggests 
that IBCPs are indeed ubiquitous in biological systems p[,|lG| . The abundance of non-bilayer lipids in biological systems 
is a longstanding puzzle, and several studies have explored the possible relationship between lipid polymorphism and 
biological membrane function (for reviews see Refs. ||ll|-|l|l). In particular, it was shown that certain cells like E. 
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Coll regulate the relative amount of bilayer and non-bilayer lipid in order to control membrane morphology and 
function of membrane proteins. E. Coli mutants lacking the non-bilayer lipid phosphatidylethanolamine were found 
to be severely impaired in regard to protein transport across the plasma membrane . Luzzati and coworkers have 
pointed out that if any IBCP had a biological function, it should be stable at physiological temperatures and in excess 
water IBCPs are also relevant for biotechnological applications, e.g. they have been used recently as artificial 

matrix which enables membrane proteins to crystallize in a three-dimensional array . 

Since IBCPs consist of a lipid bilayer separating two labyrinths of water, is it reasonable to assume that the two 
monolayers are arranged in a way which is locally symmetrical about their mid-surface. In differential geometry, a 
surface which curves to both sides in the same way is known as a minimal surface, i.e. a surface for which the mean 
curvature vanishes at each of its points. Therefore the mid-surfaces of the IBCP's lipid bilayers are often modeled 
as triply periodic minimal surfaces (TPMS, sometimes also denoted IPMS for infinite periodic minimal surfaces) of 
cubic symmetry [ p^ , pT[ . Indeed it has been confirmed by a thorough analysis of electron density maps derived from 
X-ray data that the mid-surfaces of these structures can be very close to minimal surfaces [|] . Note that the stronger 
requirement of the monolayers' global symmetry about the mid-surface leads to a certain subclass of cubic TPMS 
which is called balanced since then the two water labyrinths have to be congruent to each other. 

Since the seminal work of Schoen jl8|, it is known that there exists a large variety of different cubic TPMS. Before 
Schoen's work, only P, D and C(P) were known from the 19th century work of Schwarz and his students; then 
Schoen described G, F-RD, I-WP, 0,C-TO and C(D). More cubic TPMS have been found later by Karcher and 
Polthier Fischer and Koch and others, but most of them seem to be too complicated as to be of physical 
relevance. In Figs. |^ and ^ we depict surfaces corresponding to the seven different TPMS treated in this work. Until 
now, only G, D and P have been established to be realized as IBCP in lipid- water mixtures One system for 

which all three of these phases is stable is 2:1 lauric acid (L A) /dilauroyl phosphatidylcholine (DLPC) and water 
(with P coexisting with excess water) The only other system for which G, D and P are known to be stable is 

didodecylphosphatidylethanolamine and water [ p2| ; however, here G exists with hardly no hydration at all, and the 
regions for D and P are very difficult to access experimentally, thus this system should not be compared with the 
theoretical results for fully hydrated phases presented here. In all other systems with IBCPs, only one or two of them 
are stable (for example, for monoolein and water, G and D are stable, with D coexisting with excess water [p3[). In 
this work we present a model which can explain why only G, D and P are observed in lipid- water mixtures and why 
their actual occurrence seems to depend on the specific system under investigation. The phase diagram predicted by 
our model agrees well with the experimental phase diagram for 2:1 LA/DLPC and water. 



Exact (Weierstrass) representations are only known for G, D, P and I-WP [g4-27|. However, it has been shown 
in recent studies |^,^ that these and other cubic TPMS can be generated as isosurfaces of a density field, which 
minimizes the free-energy functional of a simple Ginzburg-Landau model for amphiphilic systems pO) . Using the 
Fourier representations of Ref. (29|, it is now possible to investigate also the physical properties of the IBCPs, for 
whose minimal mid-surfaces no exact representations are known. In particular, the question of whether these structures 
might correspond to stable IBCPs in lipid-water mixtures can now be addressed quantitatively. In this work we study 
G, D, P, I-WP, S, C(P) and F-RD, of which G, D, P, S and C(P) are balanced. 

The most important contribution to the free energy of amphiphilic interfaces is the curvature energy, which dom- 
inates all other contributions if the radii of curvature are large compared with molecular length scales |^l[. For 
the inverse phases in lipid-water mixtures, it is generally assumed that another important contribution to the free 
energy is the stretching energy of the hydrocarbon chains ||3^~|3^ . In order to attain optimal free energy, the neutral 
surfaces of the two monolayers should realize at the same time constant mean curvature cq and constant distance I 
to the minimal mid-surface. Since this is not possible geometrically, the free energy will always be frustrated 
The relative importance of interface (bending) and bulk (stretching) contributions in self-assembled interfaces is a 
long-debated issue. For amphiphilic systems, for which the free energy contributions come mainly from the interfaces, 
the main problem is of a mathematical nature: the overall solution can be investigated only in a model which allows 
the interfacial positions to adjust themselves freely. In contrast, in diblock copolymer systems bulk contributions are 
much more important and models based on random polymer coils and interacting monomer densities are more ade- 
quate. Self-consistent mean-field theory for diblock-copolymer systems suggests that in this case, interfacial tension 
and stretching are equally important |36|] . Recently, this framework has been adapted to treat lipid-water mixtures 
prf . This approach is complementary to the membrane approach; it is most useful for the investigation of systems, in 
which the internal structure of the bilayer is not uniform. However, self-consistent field theory for large molecules is 
not expected to describe amphiphilic systems well at large lattice constants, when the membranes' internal structure 
is essentially preserved and the main contributions to the free energy depends on their shape. In particular, only in 
the membrane approach does it becomes possible to address the role of thermally activated membrane undulations. 

In order to study interface-dominated systems, it has been suggested in Ref. to consider the two classes of 
surfaces, which can completely relax one of the two relevant contributions. Parallel surfaces to the minimal midsurface 
realize the optimal chain length I and are frustrated only in their bending energy, while constant-mean-curvature 
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surfaces (CMC-surfaces) realize the given spontaneous curvature, but are frustrated in their stretching energy. The 
mathematical properties of parallel surfaces are well known in terms of the properties of the underlying minimal 
surface; the parallel-surface model has therefore been the main tool for the investigation of IBCPs [ ^2|j3^j39| , |2l[ | . Up 
to now, the main difficulty of this approach has been to take into account the variation of Gaussian curvature over the 
minimal mid-surface. In earlier work, the curvature energy was therefore expanded in moments of the distribution of 
the Gaussian curvature. In Re f. ||3^ this was done up to eighth order for D, and in Ref. up to second order for 
P, D and G. CMC-surfaces 34[2^ have the advantage that they offer an unified picture for all relevant phases — 
lamellar, inverse bicontinuous cubic, inverse hexagonal and inverse micellar cubic phases are modelled by planes, cubic 
CMC-surfaces, cylinders and spheres, respectively. However, since no exact representations are known, the stretching 
frustration for IBCPs has so far been evaluated numerically for D only [Q. 

We will show below that in the case of short chains (or large lattice constants) , stretching energy is prohibitively 
large, and the two monolayers of a IBCP in lipid-water systems can be modeled by parallel surfaces. In the case of 
long chains (or small lattice constants), we expect stretching and bending to be equally important; this agrees with the 
results for diblock copolymers, which are in the long chain limit. Our reasoning implies that CMC-surfaces are not a 
good representation of IBCPs in lipid-water systems for physical reasons. However, for small distances to the minimal 
mid-surface, their mathematical properties are very similar to the ones of parallel surfaces thus many results 

from the parallel surface model will carry over to the CMC-model. In general, however, the parallel surface model is 
more appropriate both for physical and mathematical reasons. Moreover, there is strong experimental evidence for 
it, both from a detailed reconstruction of electron densities |^ and from analysis of swelling data pl|| . 

In contrast to earlier work, we evaluate the curvature energy in the parallel surface model without any expansion and 
to high numerical accuracy by first calculating the distribution of Gaussian curvature over the minimal mid-surface. 
For G, D, P and I-WP this can be done from their Weierstrass-representations. We also consider the structures S, 
C(P) and F-RD by calculating the distribution from the representations which we obtained from Ginzburg-Landau 
theory. Our main result is that the free-energy densities of the structures G, D and P are considerably lower than 
those of the other investigated structures due to their narrow distribution of Gaussian curvature over the minimal 
midsurface. We show that this result persists when thermal membrane undulations and higher order terms in the 
bending energy are considered, and argue that potential additional contributions to the free energy are unlikely to 
change it. This explains why only G, D and P have been observed in lipid-water mixtures. In fact, the calculated 
phase diagram agrees nicely with the experimental one for 2:1 LA/DLPC and water, for which these three phases 
coexist. We show that due to the existence of a Bonnet transformation between G, D and P, these phases coexist 
along a triple line in our model. Simultaneously, P coexists with an excess excess water due to a mechanism called 
emulsification failure. The Bonnet transformation also implies that the free-energy densities of G, D and P scale as 
a function of concentration with a universal geometrical quantity, which we term topology index. Since the topology 
index decreases from G to D to P, the gyroid G is most prominent, followed by smaller regions of stability for D and P 
at higher water concentrations. Any additional contribution to the free energy, which introduces a new length scale, 
is expected to change the delicate balance between these three phases; this includes stretching contributions, van der 
Waals and electrostatic interactions. Since such additional effects are specific to a given experimental system, this 
result explains qualitatively why usually only one or two of these phases are observed in lipid-water mixtures. 

The paper proceeds as follows. In Sec. |lj, we specify the free-energy expression which has to be evaluated for 
each type of minimal mid-surface, and discuss some imm ediate consequences of the model. The calculation of the 
distribution of Gaussian curvature is explained in Sec. HI for P, D, G and I-WP from their Weierstrass representation. 
Results for the Gaussian-curvature distributions for S, C(P) and F-RD are obtained in Sec. iV|from the representations 
derived recently from a simple Ginzburg-Landau model. These results are combined in Sec. ^to predict phase behavior 
and to numerically calculate phase diagrams. In Sec. Vl, we discuss which physical mechanisms will bre ak th e Bonnet 
symmetry between P, D and G. Finally, we compare our results with experiments and conclude in Sec. VII. 



II. THE MODEL 

The bending energy of one lipid monolayer is described by the Canham-H elf rich Hamiltonian 

Eb = J dA' {2k {H^ ^ co)^ + K K^} , (1) 

where the integration extends over the monolayer's neutral surface, which is at a distance I from the minimal midsurface 
of the bilayer. In the framework of the curvature model, bending and stretching contributions decouple at the neutral 
surface [ pl| ; this allows to disregard area-stretching contributions. The neutral surface has been shown experimentally 
to be located close to the polar-apolar interface Q]. In the curvature model (Q), the neutral surface is characterized 
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by mean curvature Gaussian curvature K and differential area element dA\ the superscript / indicates that the 
neutral surface is a parallel surface. The model parameters are the spontaneous curvature cq, the bending rigidity n 
and the saddle-splay modulus R of the monolayers. For lipids, the bending rigidity k for a monolayer is of the order 
10 — 20 fcsT. Usually the saddle-splay modulus R is assumed to have a small negative value; in fact the curvature 
model Eq. (|l|) without any additional constraints is well-defined only for —2k < k < The spontaneous curvature 
Co can be considered to result from the mismatch between optimal head group area, volume incompressibility of the 
chains and optimal chain length. For non-ionic lipids, it depends mainly on temperature, which changes the head 
group hydration and therefore optimal head group area. Since lipids typically form lamellar and inverse phases, where 
the monolayers prefer to bend towards their polar sides, we define positive curvature to be towards the water regions. 
We assume that spontaneous curvature scales with temperature as cq {T — Tb), since such a linear relationship 
seems to hold very well for ternary surfactant systems |Q . 

The stretching energy can be assumed to be harmonic about the average chain length I, so that 



Es= J dA ks (L - (2) 

where L is the local chain length. The stretching modulus ks is itself a function of the average chain length /. The 
relevant scaling law has to be obtained from a microscopic model; here, we model the lipid monolayer as a Gaussian 
polymer brush grafted to the polar-apolar interface, but more realistic models give similar results | |3l| . The stretching 
energy per chain reads Eg = {iksT /2N){1 / , where uk is the (microscopic) Kuhn length and N the number of 
Kuhn segments. Enforcing the volume constraint Id^ — Na\- (where is the grafting density), and transforming 
from energy per chain to energy per area, results in a potential of the form Eq. (H) with kg ^ kBT{aK /Id'^). Thus 
the stretching modulus kg scales inversely with average chain length / — stretching becomes very difficult for short 
chains |^ . 

Although the bending energy Eq. (|^) is well known not to depend on microscopic details, for our purpose it is 
important to note that the bending rigidity k also depends on the average chain length I. Starting from the stretching 
energy of Eq. (|^) and taking into account how the chain volume varies with curvature, results in the bending energy 
of Eq. (|]) with K ^ kgl* ^ kBT{aKl^ / d^) The result that the bending rigidity k scales with the third power of 
film thickness I is well known also from the elasticity of thin solid sheets. 

We now can estimate the relative importance of bending and stretching contributions. It is known that for both 
CMC-surfaces and parallel surfaces, the normalized standard deviations of L and H , respectively, are almost inde- 
pendent of hydrocarbon volume fraction 132], so that {{L — Z)^) ~ P and {{H — (H))^) ~ (H)'^ ~ a~^, where a is the 
lattice constant. Therefore Es/Eb ~ ksr/Ka~'^ ~ [a/lY. Thus, for small chain length I (large lattice constant a), 
the stretching term is much larger and the chains can be assumed to be of constant length. This corresponds to the 
well-known result that the curvature energy is the relevant contribution to the free energy as long as the curvature 
radii (that is the lattice constants) are much larger than any molecular length |3^]. Only for large chain length I 
(small lattice constant a), the bending energy becomes comparable and both terms have to be considered. Experi- 
mentally, this corresponds to the fact that for small water content (i.e. chain length I of comparable size as lattice 
constant a), the bicontinuous cubic phases become unstable with respect to the hexagonal phase, whose geometrical 
properties require larger chain stretching. Note that the relative importance of stretching and bending is asymmetric: 
since the chain length I cannot be larger than the lattice constant a, there is no regime in which bending can be 
neglected. Thus CMC-surface will never be a good approximation for monolayers in lipid-water systems. However, 
the stretching frustration can be relaxed by swelling the bilayer with oil; therefore CMC-surfaces are good approxi- 
mations for amphiphilic monolayers in ternary systems with oil, water and amphiphile po|] . Similar considerations of 
the relative importance of stretching and bending can be found also in the context of large membrane inclusions (like 
transmembrane proteins) where chains have to stretch in order to decrease the hydrophobic mismatch at the inclusion 
boundary Q. 

For the rest of this paper, we proceed within the framework of the parallel surface model. We assume that the 
interfaces between the tail and head group regions of the two lipid monolayers in a ICBP are located at distances ±Z 
away from a given cubic TPMS; typically I has values between one and two nanometers. Note that for a TPMS the 
mean curvature satisfies H — (ci -|-C2)/2 — everywhere (where ci and C2 are the two principal curvatures), while the 
Gaussian curvature K = ciC2 = — ci^ varies as function of the position and is only restricted to satisfy K <Q. The 
differential area element dA\ the mean curvature and the Gaussian curvature on the parallel surfaces then 
follow as functions of I and the quantities dA and K on the minimal mid-surface | p^ , |3^ , with 

dA^^dAH^KP), i^'^^, ^'-YVm^- (3) 

Note that these formulae are special cases of Steiner's theorem of integral geometry; they are therefore no approxima- 
tions for small I, but exact. Since minimal surfaces have K <Q, positive mean curvature is defined here to correspond 
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to positive /. 

The effective bending energy of the hpid bilayer is only a function of the Gaussian curvature K of the minimal 
midsurface since its mean curvature H vanishes. In an expansion for small chain length it can be expressed as an 
infinite series in powers of K. Using Eq. (0) and Eq. (||) yields up to order P 

Ef,^ jdA{AclfK+{2K + %CQlK + AclfK)K + AKpK^] . (4) 

Therefore the effective saddle-splay modulus for the lipid bilayer, = 2k + ScqIk + Ac^Pk, is not just the sum of 
the monolayer's moduli, but is corrected to higher positive values due to the presence of the spontaneous curvature 
Co 0- We conclude that as long as cqI > — k/4k, the IBCPs are favored over the lamellar phase, since the preferred 
curvature of the monolayers translates into a topological advantage of saddle-type bilayer structures. The third term 
in Eq. (^) favors saddle-type structures in any case, and also removes the degeneracy which minimal surfaces might 
experience under deformations which preserve the minimal surface property. This has been noted in a seminal work 
by Bruinsma |^^; in his notation, we have R = AkP for the elastic modulus of the K^-teim. 
The volume fraction occupied by the hydrocarbon can be calculated as 

where we have used the Gauss-Bonnet theorem J dAK — 2Trx- A* and x denote surface area and Euler characteristic 
of the given cubic TPMS in the conventional unit cell with unit lattice constant. For the following, it is useful to 
introduce the topology index T — {A*^ /2tt\x\Y^'^ ■ This quantity is independent of scaling and choice of unit cell 
and characterizes a given TPMS's topology in an universal way; the higher its value, the smaller is the porosity 
and the larger the specific surface area of a TPMS. It is the only quantity which characterizes a two-dimensional 
minimal surface in three-dimensional space independently of lattice constant and choice of unit cell, and its relevance 
for the description of bicontinuous cubic phases in amphiphilic systems has been discussed before [ ^ , ^9|j2g|j40| , ^ . 
The topology index of the various TPMS is of order 1 and is given in Tab. |. The gyroid G is expected to have the 
highest value since it divides space into labyrinths with threefold coordinated vertices. S has a similar value for F, 
since it features a mixture of three- and fourfold coordinated vertices. Then comes the double diamond structure 
with fourfold coordinated vertices, while all other structures have higher values. In Tab. | we collect the values of 
A*, X and F for the cubic TPMS investigated in this work as obtained from our numerical representations in the 
framework of a siniple Ginzburg-Landau theory (see below). Exact values are known for G, D, I-WP, P and C(P) 
and given in Tab. ^ Comparing numerical and exact values shows that for most structures, our numerical results 
are quite accurate. Due to numerical limitations for the more complicated structures, the values for S and F-RD are 
less precise. In particular, we expect S to have a smaller value for F than G. 

Equation (^) can be inverted numerically to give the dimensionlcss function a/ 1, the lattice constant a in units of 
the chain length Z, as a function of hydrocarbon volume v. For small v we find 

In Fig. |l| we plot the first-order approximation a/l = 2A* /v of Eq. (^) for all structures considered; it can be hardly 
distinguished from the full curve over the full range of v. By visual inspection of our numerical representations (see 
below), we find that for all structures considered, parallel surfaces do not self-intersect as long as w ^ 0.8; IBCPs 
can accommodate large amounts of hydrocarbon since their geometry is so close to that of the lamellar phase. A 
typical value for I is ISA, so that for v = 0.5 we have as = 325 A, uf-rd — 285 A, = 230 A, acp = 211 A, 
ai-wp — 208 A, ao = 185 A and ap = 141 A. Note that at a given v, the ratio of the lattice constants of two 
coexisting IBCPs is simply the ratio between their scaled surface areas A* . In Figs. || and || we depict for each of 
the structures considered one of the two monolayers for the hydrocarbon volume fraction v = 0.5. Then the volume 
fraction of each of the two labyrinths is (1 — v)/2 = 0.25. 

We are now in a situation to rewrite the curvature energy from Eq. (^ for our purpose. We consider a typical experi- 
mental situation, in which temperature is controlled at constant volume. As explained above, changing temperature T 
amounts to changing spontaneous curvature cq. Under the assumption of incompressibility of both components, there 
is only one degree of freedom for composition, which we take to be the hydrocarbon volume v (which is equivalent to 
the amount of interfacial area). Dividing the free energy by the constant volume, we arrive at the free energy per unit 
volume as a function of v. However, v is not a controlled quantity, since at two phase coexistence, it is not clear how 
it will distribute into the two phases. The controlled quantity can be considered to be the chemical potential coupled 
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to V (that is the derivative of free energy - or free energy volume density - for hydrocarbon volume fraction - or 
interfacial area), and therefore phase transitions have to be calculated by using the Maxwell construction for the free 
energy density as a function of v (which amounts to requiring the same chemical potential for both phases) . We want 
to remark parenthetically that the free energy per unit area cannot be used to calculate phase coexistences, since an 
exchange of area, controlled by the chemical potential of the lipid, has to be allowed between coexisting phases. 

We write the free energy volume density dimensionless by using a factor 1/{2kcq). For the lamellar phase the 
hydrocarbon volume fraction v — 2l/a, thus its dimensionless free-energy density follows from Eq. (|l|) to be 

The larger v, the smaller a and the more frustrated curvature energy per volume accumulates. 

For IBCPs, we again use Eq. (|l|) and Eq. (^ in order to write the dimensionless free-energy density as a function 
of hydrocarbon volume fraction v and the properties of the bilayers' minimal mid-surface: 




(8) 

where again we have used the Gauss-Bonnet theorem and defined r — —k/2k. As mentioned above, the curvature 
model Eq. (^ without any additional constraints is well-defined for < r < 1 and bicontinuous phases are favored for 
cqI ^ r/2. A*, K* and x refer to a conventional unit cell. Using Eq. (^) to lowest order, a/l — 2A*/v, and defining 
E{K*) = K*A*/8ttx, we rewrite Eq. (|) as 

dA* /«\2\"V, l+fo„. 

Co 



In Eq. (^) and below, spontaneous curvature cq is measured in units of l/l. Since for the minimal midsurface both 
K* and x are negative, 'E.{K*) is a positive quantity which varies over the surface. Note that J {dA* / A*)'E.(K*) = 1/4 
due to the Gauss-Bonnet theorem. The free-energy density Jl of the lamellar phase, Eq. (|^), is a special case of /(,, 
with 'B.{K*) = X = 0- In this work, we will treat the free-energy density Eq. without any further approximations. 
Previously, similar expressions have been expanded in small v (which is equivalent to small I) ]3^,|39| , similarly as 
discussed in Eq. (^). In fact, each term (KP)" in an expansion for small I corresponds to a term v^nTi jj-^ g^j^ expansion 
for small hydrocarbon volume v. Since KP = —E{K*){v/r)^ , higher-order terms in Eq. (|l|) would result in a free- 
energy density which has a similar structure as Eq. (^) . As we will see below, this means that our results will remain 
valid in the case that higher-order terms are included in the curvature energy of the monolayers. 

It is well known that thermal fluctuations can contribute significantly to the free-energy of amphiphilic systems. 
Since we consider the regime where the stretching energy dominates, undulations are the favored fluctuation modes. 
For the lamellar phase, they give rise to steric repulsion EG], which in our case can be written as 



Cn 



knT 



f steric ~ oo 2 ( ) /i \? (^0) 

32Cq \ K J (1 — w)^ 



where Coo = 0.106 from Monte Carlo simulations ||51|,|5^ and field-theoretical calculations ||5^]. For an IBCP, steric 
repulsion is hardly relevant since the lateral correlation length ^|| is of the order of the lattice constant a. Therefore, 
the perpendicular correlation length ~ y^ksT/K ^|| ^ ^|| ^ a, and few membrane collisions should occur. However, 
fluctuations of the lipid bilayer lead to a renormalization of n and R, for all phases. Since k multiplies the average 
mean curvature squared, which vanishes for the mid-plane, only the renormalization of k at length scale i has to be 
taken into account. We identify the typical length scale of a cubic structure with (K)^^/^. This implies i/l = 2r/w, 
so that r in Eq. (p|) gets renormalized, with 



rR 



5 fcsT, f2T 
■ In 



127r K \ V 



(11) 



Since the renormalization increases R, this effect increases the topological advantage of the saddle-shaped structures as 
does the spontaneous curvature. Note that for the lamellar phase, these effects might favor the formation of wormholes 
psf . Here, we neglect this aspect, as well as the acoustic modes of the bicontinuous structures at large wavelengths, 
which have been discussed by Bruinsma E^. Note also that the thermal fluctuations also give a contribution to the 
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free energy, which is hnear in the membrane area. We have omitted this term here, since it can be absorbed into the 
chemical potential of the lipid. 

Having defined the model, we now turn to its phase behavior. We first note that the dimensionless mean curvature 
averaged over the parallel surface follows from Eq. (^) as 

where again l/a = v/2A* and the Gauss-Bonnet theorem has been used. It is plotted in Fig. ||a. Note that all curves 
fall on a universal curve when the hydrocarbon volume fraction v is scaled with the topology index F. It reaches the 
value Co at 



4co 
1 + co 



(13) 



Thus, Co is realized in a sequence given by increasing topology index; in particular, the gyroid structure G reaches Cq 
at the highest value of v. We conclude that if all considered structures were stable, they would appear in the sequence 
G - S - D- I-WP - P etc with increasing water concentration. 

The stability of the different phases is determined by the free-energy density given in Eq. (||). The second (topologi- 
cal) term is easy to understand. The topological properties of the parallel and mid-surfaces are the same: for negative 
saddle-splay modulus k (r > 0), the different IBCPs are favored according to their values of the topology index F. 
Here the G-structure performs best since it has the lowest porosity. The first term in Eq. (||) is more complicated, 
since the mean-curvature properties of the parallel surfaces translate into the Gaussian-curvature properties of the 
minimal mid-surface in a quite complicated way. As discussed above, the spontaneous curvature of the monolayers 
implies a topological advantage of the bicontinuous phases. For a more detailed analysis, which includes all terms of 
the expansion, we note that the first term in Eq. (H) measures the standard deviation of the parallel surface's mean 
curvature from the spontaneous curvature. Since it follows from Eq. (|^) that the distribution of mean curvature 
over the parallel surface, iJ'dA', is proportional to the distribution of Gaussian curvature over the minimal surface, 
KdA, we expect that those structures will be more favorable which have small standard deviations for their Gaussian 
curvature distributions p9[ . In order to quantify this concept, we now turn to the distributions of Gaussian curvature 
which later will allow us to evaluate Eq. (|§) without any further approximations. 



III. CURVATURE PROPERTIES FROM WEIERSTRASS REPRESENTATIONS 



Weierstrass representations are known for P, D and G |25y27]] and I-WP ||24| , |26| . For each of these TPMS, a 
fundamental domain can be identified, so that the rest of the surface follows by replicating it with the appropriate 
space group symmetries (/m3m, Pn3m, laZd and ImSm, respectively). The Weierstrass representation is a conformal 
mapping of certain complicated regions within the complex plane onto the fundamental domain: 



{xi,X2,X3) — Re / dz 
Jo 



Riz) {1- z',i{l + z'),2z) (14) 



where (m, v) arc the internal (and conformal) coordinates of the minimal surface. The geometrical properties of such 
a surface follow as 

dA{z)^\R{z)\^ {l + \zff dudv, H{z)^0, K{z)^ (15) 

l-R(^)l (1 + kl ) 

with z = u + iv. Obviously the (isolated) poles of R{z) correspond to the flat points {K = 0) of the minimal surface. 
Only few choices of R{z) yield embedded minimal surfaces. The ones for D and P have been known since the 19th 
century from the work of Schwarz: for D it is R{z) = {z^ — 14z^+ 1)^2 . p corresponds to the same region and follows 
simply by the Bonnet transformation R{z) e^^R{z) with 9 — 90°. Equation ( [l5| ) implies that P and D have the 
same metric and the same distribution of Gaussian curvature. However, since they map differently into embedding 
space, they have different space groups and lattice constants. The gyroid G was discovered in 1970 by Schocn |lj] 
as another Bonnet transformation of D, with 6 = 38.015°. The Weierstrass representation for I-WP was found only 
recently |2j,|2^. If one of its poles is chosen to be at infinity, one has R{z) = {z{z^ + l))-i. In Fig. | we show 
the fundamental domains in the complex plane of the different Weierstrass representations. In Tab. || we collect the 
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values for scaled surface area A* and lattice constant a as they follow from the Weierstrass representations. We also 
give the number N of replications of the fundamental domain needed to build up the surface in one conventional unit 
cell and the values for the Euler characteristic x which follows from the topology generated in this process. Note that 
for D one can choose another conventional unit cell which is contained n ~ 8 times in the one chosen here. This has 
been done in Rcfs. then %, A and K scale as n~^, n~^l^ and n^^/"^, respectively. In particular, one then has 

X = —2 and A = 1.9188925. Our choice is motivated by the fact that it is the relevant one also for ternary systems. 
We define the distribution function of Gaussian curvature, 

f{K) = J dA{u, v) 5{K ~ K{u, v)) , (16) 

over some surface parametrized by internal coordinates {u, v), where in our case dA(u, v) and K(u, v) are determined by 
Eq. (p^. For the generating functions R{u, v) detailed above, it is not possible to calculate / analytically. Therefore, 
we evaluate it numerically by considering a histogram {fi} for a set of discrete Gaussian curvature values {Ki} spaced 
equidistantly with AK. We cover the (m, u)-region corresponding to the fundamental domain with a square grid of 
M X M points. For each point, we calculate dA and K from Eq. (jlj) and add dA to the fi which corresponds to 
the Ki with K G [Ki — AK/2, Ki + AK/2]. In order to obtain results for one conventional unit cell, the values {fi} 
are then multiplied by the number N of replications and the values {Ki} are scaled with (both N and a are given 
in Tab. ||). In practice we use M = 30000 and AK = 0.025; the latter value amounts to about 150 bins for the 
histogram. By this procedure, we have 

/, = / dA{u, v) / dK' 5{K' ~ K{u, v)) = AKf{K,) . (17) 

J JKi-AK/2 

The moments Sn of the distribution f{K) follow as 



Sn — 



J dKf{K)K^ = J dAK^ - f^^^"" ■ (18) 



In particular, Sq Y.i fi = ^* a-i^d ■^'i = Y.i h^i = 2ttx- 

The procedure described here was carried out to obtain {Ki, fi} for D and I-WP. Since P and G are related to D by 
a Bonnet transformation, their distributions can be obtained simply by appropriately rescaling. For example, for the 
transformation D ^ G, the Ki have to be rescaled with a^/a\, and the fi with A*q/A*j-,, as given in Tab. ||. In Fig. |a, 
we plot the four different distributions obtained. P, D and G feature one peak around a (relatively large) K-value, 
that means they are quite uniformly curved, although they have flat regions, too (which correspond to the poles of 
R{z)). I-WP is bimodal, that is it has prominently both flat and strongly curved regions. The strong curvature of 
D is somehow artificial due to our choice of unit cell. If one corrects for this by a factor 1/4 for the K- values, it 



becomes clear that I-WP has a much broader distribution than P, D and G. This can also be inferred from Tab. II], 
where we list the first moments Sn of the distributions for D and I-WP. The ones for D are identical to the ones 
given in Ref. after appropriate rescaling Sn with 8"*^^^+^'/'^ due to the different choice in unit cell (the numerical 
deviation is less than 0.02 % even for TV = 8). The moments for P and G can be calculated from the ones for D due to 
the Bonnet-transformation between them; e.g. S'^ = S-^{Nc/N]j){ao/aD)'^'"'^~^'' with the values for TV and a given 
in Tab. In particular. A}. = A*jj{NGa]j) / {Noal), xg = Xd{Ng/Nd) and Pg = TD{NGal)/{NDal) — as can be 
verified from Tab. It follows immediately that the quantity ^(K*) — K* A* /8ttx defined for Eq. (g) is invariant 
under a Bonnet transformation. 

In order to quantify the variance of the different distributions, we calculate standard deviation, A — {{K — 



{K)Y)/{KY = S0S2/SI - 1 (here (...) means area average). This 



cell and Bonnet-transformation, and given in the last row of Tab. Ill 
for P/D/G. 



uantity is independent of scaling, choice of unit 
. As expected, its value for I-WP is higher than 



IV. CURVATURE PROPERTIES FROM GINZBURG-LANDAU 



In order to investigate the curvature properties of S, C(P) and F-RD, we use the representations as isosurfaces to 
fields <I'(r) which we obtained recently as local minima of a simple Ginzburg-Landau theory for ternary amphiphilic 
systems |29| . The same model can also be employed for binary systems if the scalar order parameter field <I>(r) 
is interpreted to be the local concentration difference between water in the first labyrinth and water in the second 
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labyrinth. Then the mid-surface of the hpid bilayer is identified with the surface ^{r) = 0. The free-energy functional 
is given by 

J^m^ J dr {(A<i>)2+5(ci>)(V<i>)2 + /($)} . (19) 

A reasonable choice for / and g has been found to be 

/($) = + + j^) ^ = 50 + .g2$2 . (20) 

For given model parameters (.goi52, fo), the free-energy functional is minimized in Fourier space by implementing the 
correct space group symmetry and minimizing for Fourier amplitudes and lattice constant with conjugate gradients. 
For balanced TPMS, an additional black and white symmetry has to be implemented. The resulting fields <i>(r) 
correspond to the various local minima. We showed in Ref. that by tuning (go, g2, /o) appropriately, the <i>(r) = 0- 
isosurfaces can be made to become very close to minimal surfaces, and generated representations for nine different 
TPMS: G, D, P, S, C(P), C(D), C(Y), I-WP and F-RD. 

K-distribution can be obtained from these representations as follows. By triangulating the isosurfaces with the 
marching cube algorithm, we numerically obtain a metric on the surface. From the Fourier representations of the 
fields $(r), it is possible to calculate exactly the K- values on the vertices of the triangulation. Using a histogram 
technique like in the Weierstrass case, it is thus possible to obtain {Ki, fi} even for structures like S, C(P) and F-RD, 
which might be of physical relevance but for which no exact representations are known. The distributions obtained in 
this way for go — —3, g2 — 7.01 and fo — are plotted in Fig. |5|3. All three of them are roughly trimodal and feature 
much stronger curved parts than the ones discussed above. This reflects their complicated shape, which makes them 
more difficult to access numerically than the ones which have been generated by Weierstrass representations; it also 



explains the difficulty to obtain smooth data for the K-distributions. In Tab. Ill, we give their first moments Sn 
and the measure A of their variance. All of them have broader distributions than I-WP, in the sequence S - F-RD - 
C(P). Therefore we can conclude already from the distributions that the IBCPs should become less favorable in the 
sequence P/D/G - I-WP - S - F-RD - C(P). 



V. PHASE BEHAVIOR 



Given the distribution of Gaussian curvature as histogram {Ki,fi}, it is straightforward to evaluate the first term 
of the curvature energy, Eq. (^), of the bicontinuous phases. Since it is an area average over a complicated function of 
K* , we have to replace J dA* and K* by ^ fi and Ki, respectively. In Fig. we plot the quantity {{H^ ~ cq)^)iP as 
a function of hydrocarbon volume fraction v for cq = 1/6. In contrast to {H^)il, compare Eq. (^2|), no simple formula 
exists, and it has to be calculated numerically as 



+ Co 



2 



((^ - ^'or)ii' = — - — ^—^ — ■ (21) 



Since is invariant under a Bonnet transformation, {{H^ — co)^);P is the same for P, D and G, except for a rescaling 
of V with r. In particular, the curves reach the same minimal values and the minima are shifted to the right with 
increasing topology index F. Thus the structures P, D and G perform equally well, but at different values of v. From 
Fig. raa, we see that all other structures cannot reach this low value of frustration since they have much broader 



distributions than P, D and G (compare Tab. ill). Note that the scaling with v/T is valid in Eq. ( |12| ) for all structures 
and in Eq. ( ^l|) only for those related by a Bonnet transformation. 

Except for a factor of v/cq, the free-energy density (^ for r = is identical with ((iJ' — coY)iP (see Eq. (pi])). 
In Fig. 1^, we show the free-energy density of all considered phases for cq — 1/6, r = and n/ksT = 10. ^From 
the graph, we infer the phase sequence L - G - D - P with increasing water content; the other phases are not stable 
since they cannot reach that small an amount of bending frustration as can P, D and G. As discussed above, Eq. ( [T^ ) 
implies that the sequence G - D - P is determined by the topology index. 

Since the free-energy density for P rises again for small w, it would be more favorable for it to coexist with excess 
water. This mechanism is known as emulsification failure in surfactant systems and allows to extend the sequence 
to include a P-l-W coexistence. The emulsification failure can be interpreted as a Maxwell construction of P with the 
pure water phase with / = at w = 0. It follows from Eq. (^ and the scaling of ((iJ' — cq)^);/^ with v/T that the 
free-energy density for P, D and G can be written in the scaling form 

fb = vn{v/r) (22) 
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with a universal function il{x). For small r, this function has rj(O) = 1, a minimum at a; ~ 4co/(l + cq), and diverges 
at a; ~ 4. A Maxwell construction shows that these structures always lie on a triple line, irrespective of the values for 
Co, r and k (compare inset of Fig. |^). This means that if we calculated the whole phase diagram with the Maxwell 
construction, D and P would be stable only along a line. In the following we identify phase transitions between 
G, D and P with the intersections of the free-energy curves, since under experimental conditions, small additional 
contributions to the free energy are certain to destroy the delicate free-energy balance resulting from the Bonnet 
transformation, and then will lead to extended one-phase regions for D and P (or remove one or both of them from 
the overall phase diagram) . 

In Fig. ^ we use this reasoning to construct phase diagrams as a function of v and cq for different values of r and 
K/ksT. With increasing water concentration pw = 1 — f we always obtain the sequence L-G-D-P - W-l-P. 
The most stable IBCP is the gyroid structure G close to the lamellar phase at high values of v since it has both a 
narrow distribution of Gaussian curvature and the highest topology index. At lower concentrations of v, one finds 
narrow regions of stability for D and P since they have narrow distributions, too, but increasingly lower values of F. 
Increasing cq shifts the loci of least frustration to higher v as evident from Fig. Decreasing R (i.e. increasing r) 
disfavors the IBCPs and increases the stability region for the lamellar phase, as can be seen by comparing Fig. ^ 
and ^D, where r is increased from 0.1 to 0.5, respectively. However, for cq ^ r/2 there always will be a region of 
stability for the IBCPs, even for r = 1. This stands in marked contrast to similar modelling for ternary systems where 
bicontinuous cubic phases are expected to disappear for r > 0.55 | |40[ |. The explanation lies in the fact that changing 
the spontaneous curvature Cq of the monolayers amounts to changing the effective saddle-splay modulus of the bilayer, 
compare Eq. (Q). Thus, increasing cq can counteract the effect of making R more negative. The bicontinuous phases 
are also favored by lowering k without changing r, since this increases the fluctuation effects which favor IBCPs 
(positive renormalization of R) and disfavor the lamellar phase (increase of steric repulsion). This can be seen by 
comparing Fig. and ||c, where K/fc^T is decreased from 10 to 2, respectively, while r is kept constant at 0.5. 



VI. BREAKING THE BONNET SYMMETRY 

It is clear from our discussion in Sec. ^ that all contributions to the free energy, which can be written in the scaling 
form vn{2l{Ky/^) = vn{v/r), wiU preserve the triple-line coexistence of G, D and P. Obviously, an interaction 
which introduces an additional length scale, in addition to I and (_?C)~^/^, is needed to break the Bonnet symmetry 
and to lift the triple-line coexistence. Here we want to discuss van der Waals interactions and chain stretching. 

The non-retarded van der Waals interaction at distance d — |ri — has the general form 

w{d) = -^pip2 d-"" (23) 

where Ah is the Hamacker constant, and pi and p2 are the number densities of atoms in the two interacting bodies. 
The total interaction energy between two colloidal particles of volumes Vi and V2 is then given by 

d\i f dVa w(|ri -rsl) . (24) 

JVi JV2 

Unfortunately, the energy density of the van der Waals interaction cannot be calculated easily in a cubic bicontinuous 
phase. However, it is clear from Eq. (|2j) that it is the embedding of the TPMS in three-dimensional space, which is 
important for the van der Waals energy. Distances between different parts of the bilayer are not determined uniquely 
by the iiT-distribution, which describes the internal geometry of a surface. This is evidenced by the self-intersections 
of the parallel surfaces for v ^ 0.8, which do not affect any of the considerations of the curvature energy described in 
the previous sections, but would give a divergent contribution to the van der Waals energy. Thus, we conclude that 
the van der Waals interaction must break the Bonnet symmetry. We illustrate this point in Appendix where we 
approximate the TPMS by a FCC array of spherical shells. 

In order to estimate the frustration energy due to chain stretching, it has been suggested in Ref. to calculate 
the variance of the distance between constant-mean-curvature companions of a given TPMS. For the D-surfacc, the 
variance — {{L — Z)^) has been shown to be almost independent of v, with 

Al = aov'^a^ = 4:aDA*^f (25) 

and — 0.00035, where Eq. (^ has been employed in the second equality. Similar results can be expected for G and 
P, with prefactors ac and ap, respectively. We introduce a function a{A*), which is defined such that a{A*p) ~ ap, 
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a{A'^) = an and a{A^) = ac for the three values of A* for the P, D and G phases. Then, the stretching-energy 
density has the scahng form 



where Eg is the stretching energy from Eq. (|^) and Eq. ( p5| ) has been used for A^. Note that this energy density 
has been made dimensionless with the same factor as the curvature energy in Sec. H i.e. it is normahzed to the 
free-energy density of the lamellar phase [compare Eq^fra)]. 

These results make the scaling argument of Sec. ||inore precise. Equation ( p5| ) shows very clearly that the 
lattice constant a appears as new length scale, and that the Bonnet symmetry is broken in general by the stretching 
contributions. It is interesting to consider the possible functional dependencies of oi{A*) in more detail. First, 
if a{A*)A*'^ is a monotonically increasing (decreasing) function of A* , then P is stabilized (destabilized) by chain 
stretching due to its small value of A*, while D is disfavored (favored). In the case of increasing a{A*)A*'^ , D would not 
appear in the phase diagram, while in the case of decreasing a{A*)A*'^ , P cannot be present. Second, if a{A*)A*'^ = 
const, the Bonnet symmetry is restored. Finally, a{A*)A*'^ can of course be a non-monotonic function of A* , which 
stabilizes the phase with the smallest value of a{A*)A*'^ . 

In any case, the difference between the three phases can be expected to be small, since the prefactor ao (and 
presumably also and ap) is very small; nevertheless, the stretching energy might account for small differences, 
like the ones needed to produce the phase diagram of 2:1 LA/DLPC and water, which shows no triple-line symmetry, 
but something close to it. A more precise statement requires the calculation of for G and P in the framework 
of the CMC-model. Other estimates of the stretching energy have been discussed in Ref. they also indicate a 
breaking of the Bonnet symmetry. 



VII. CONCLUDING REMARKS 



In this work, we have investigated inverse bicontinuous cubic phases (IBCPs) in lipid-water mixtures. We system- 
atically included all IBCPs which might be expected to be of physical relevance, and found that only the ones with a 
narrow distribution of Gaussian curvature over the minimal mid-surface are stable, i.e. P, D and G. While this result 
has been anticipated in Ref. we were able to prove it quantitatively for the first time by calculating the distri- 
butions of Gaussian curvature over the minimal mid-surfaces of all relevant structures; from these we calculated the 
variance A of the ^-distributions given in Tab. Ill which is the same for P, D and G — due to Bonnet transformations 

— and higher for all other structures. 

The Bonnet transformation also implies that the free-energy densities for P, D and G have the scaling form fb = 
vD.{v/T) with hydrocarbon volume fraction v and topology index P. The topology index is a dimensionless inverse 
Euler characteristic which is independent of scaling and choice of unit cell; the higher its value, the less porous 
the structure. Since it is highest for the gyroid structure G, it has the largest region of stability. D and P have 
progressively lower values for P, and have therefore smaller regions of stability at lower values for v. 

We also used the concept of the emulsification failure (which is well established in surfactant science) to explain the 
stability of IBCPs in excess water. A Maxwell construction with the excess phase shows that P, D, G and excess water 
lie on a triple line for all values of cq and r as a consequence of the Bonnet transformation. This 'Bonnet symmetry' 
should not be interpreted as a degeneracy between G, D and P, since their different topologies ensures that their 
regions of stability are arranged according to their respective values of the topology index. It should be noted that 
the free energy per unit area (rather than per unit volume) has the scaling form f2(D/r); this might be interpreted as 
some kind of degeneracy (compare Refs. p9|pl[ |), since now all three structures attain the same minimal values. 

Our results agree nicely with experimental observations. At fixed temperature, it is always the sequence G - D 

- P which is observed with increasing water content, with G having the largest region of stability and only D or P 
coexisting with excess water ]5^ . For example, for monoolein and water at 40*^(7, with increasing water concentration 
the sequence L - G - D has been found, where G has the largest region of stability and D coexist with excess water 
psf . Our phase diagram for r = 0.1 and n/ksT — 10 (Fig. ^a) agrees particularly well with the one for 2:1 LA/DLPC 
and water [pT| (compare Fig. ^). We do not treat the inverse hexagonal phase, which experimentally occurs at large 
V and large cq, but we showed above that IBCPs will occur at large cq and small v, and argued that at large v 
stretching energy should become comparable to bending energy. This, in fact, favors the inverse hexagonal phase, 
whose geometry requires chain stretching. 

Finally, we can discuss the role of temperature for 2:1 LA/DLPC and water by using the ansatz cq — h{T — T^), 
with a balanced temperature Ti, and a phenomenological constant b. Using I — 11.4 A, we find b « 0.0011/ AK 
and Tb ss 0°C. Since the chain melting temperature is known to be T„i = 30°C, we conclude that the balanced 
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temperature Tb is well below the main transition. This seems to be the generic case for lipid systems; it mig ht be 
related, e.g., to the anomalous swelling of PC-bilayers when the main transition is approached from above |57] ]. It is 
interesting to note that the constant b obtained here for 2:1 LA/DLPC and water is very close to the corresponding 
value b = 0.0012/ AK for monolayer of C12-B5 in ternary microemulsions [^ . 

In this paper, we only considered inverse bicontinuous phases with cubic symmetry. A similar analysis for non-cubic 
phases should be straight-forward, since here Weierstrass representations are known for the main structures, especially 
for the rhombohedral and tetragonal variants of G, D and P |2^. However, we expect that these phases will not be 
more favorable than G, D and P for the same reason as for the other cubic phases. In order to make further progress 
in understanding their phase behavior, it seems more important to investigate the different mechanisms which might 
break the Bonnet symmetry between G, D and P. In our view, two calculations are needed now: (i) of the van der 
Waals energy of G, D and P and (ii) of the variance of the distance distribution of a constant-mean-curvature model 
of G and P. 



APPENDIX A: ESTIMATE OF VAN DER WAALS ENERGY DENSITY 



For a rough estimate of the energy density due to the van der Waals interaction in cubic bicontinuous phase, we 
use a FCC array of spherical shells of radius R, thickness 21 and lattice constant a. Here, the lattice constant a 
corresponds to the conventional unit cell and is taken to be the same as for the corresponding bicontinuous phase, and 
R is determined by the requirement of equal area A per unit cell, so that A* — 167r(i?/a)^. Thus we have a/l — 2A* /v 
again. The distance of closest approach between two neighboring spheres is I? = ^/2a/2 — 2R; in order to avoid 
contact between them, I? > is necessary, which implies A* < 2tt; this inequality is satisfied for all TPMS studied 
here, compare Tab. |. For the crystal, we consider the interaction between nearest-neighbor spheres only. The van 
der Waals energy between two hollow spheres is described by the Girifalco potential, which is often used to model 
the van der Waals interaction between two buckyballs Ceo [ p8| . For our purpose it is sufficient to use the Derjaguin 
approximation for the Girifalco potential, which describes the leading contributions for close approach correctly, and 
shows that the van der Waals attraction for two hollow spheres scales as ~ Ah Rl^ / . We thus find for the van 
der Waals energy density 

This energy density has again been made dimensionless with the same factor as the curvature energy in Sec. H Its 
scaling shows that two ratios of length scales, R/a ~ (A*)^/^ and a/l ~ A* /v, determine the van der Waals energy. 
It is intersting to note that the curvature-energy density in this model of spheres scales as 



co{R/l) 



(A2) 



Since R/l ^ {A*)^/'^/v ~ T/v (for fixed x), Eq. ( |A^ is consistent with the scaling function ([2^). This indicates that 
our identification of the sphere radius and lattice constant gives a consistent approximation for both van der Waals 
and bending energies. 
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0.2 0.4 0.6 0.8 y 1 

FIG. 1. Lattice constant a in units of the chain length I as a function of hydrocarbon volume v. We use a/l = v/2A* , which 
is an excellent approximation up to v « 0.8, beyond which the inverse bicontinuous cubic phases as modeled here cannot exist 
anymore due to self-intersection. From left to right, the curves correspond to P, G, I-WP, C(P), D, F-RD and S, respectively 
(the curves for I-WP and C(P) nearly collapse and therefore appear as an apparent bold line). 
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FIG. 2. One of the two monolayers of P, D, G and I-WP for hydrocarbon volume fraction v = 0.5. These structures have 
space groups PmSm, Fd3m, /4i32 and ImSm, respectively. In the full structures, a second monolayer exists parallel to the 
one shown. Each monolayer defines a space-percolating network filled with water. For the balanced cases P, D and G, the 
second network is congruent to the first; this additional symmetry changes the structures' space groups to ImZm, Pn3m and 
laSd, respectively. For I-WP, the two networks are different and only the monolayer corresponding to the I-network is shown. 
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cm 

FIG. 3. One of the two monolayers of S, F-RD and C(P) for a hydrocarbon volume fraction v — 0.5. These structures 
have space groups /43d, Fm3m and Pm3m, respectively. In the full structures, a second monolayer exists parallel to the one 
shown. For the balanced cases S and C(P), the second network is congruent to the first; this additional symmetry changes the 
structures' space groups to ImSm and laSd, respectively. For F-RD, the two networks are different and only the monolayer 
corresponding to the F-network is shown. 



16 



(a) (b) 

FIG. 4. Fundamental domains for Weierstrass representations: (a) G, D and P and (b) I-WP. Filled circles mark the poles 
of the generating functions R{z) (which correspond to the fiat points of the surfaces) and the hatched regions correspond to 
the fundamental domains of the TPMS. 
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FIG. 5. Distributions of Gaussian curvature for (a) I-WP, D, G and P as obtained from the Weierstrass representations and 
(b) S, C(P) and F-RD as obtained from our numerical representations as isosurfaces to Ginzburg-Landau local minima. Since 
P, D and G are related by a Bonnet-transformation, they differ only by simple rescaling. This data has also been shown in 
Ref. H. 
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0.2 0.4 0.6 0.8 y 1 

FIG. 6. (a) {H^)il = (i;/r)^/(4 — (w/F)^) as a function of hydrocarbon volume fraction v. The straight line is cal — 1/6. 
From right to left, the IBCPs appear in the sequence S, G, D, I-WP, P, C(P), and F-RD. (b) ((H' - cof)if as a function of 
hydrocarbon volume fraction v with co? = 1/6. The IBCPs can be identified from their sequence at the top of the figure, where 
G, D, P, S, I-WP, C(P) and F-RD appear from right to left. The straight line is {colf = (1/6)^ Since P, D and G are related 
by a Bonnet transformation, their curves differ from each other only by a simple rescaling of v with F. Therefore they reach 
the same minimal values, but at different v according to their F-values. 
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FIG. 7. Free energy densities as a function of hydrocarbon volume fraction v for col = 1/6, r = and n/kBT — 10. The 
solid line on the right is /l, the other solid lines are the different ft,. The IBCPs can be identified from their sequence at the 
top of the figure, where G, D, P, S, I-WP, C(P) and F-RD appear from right to left. With decreasing v, the phases La, G, D 
and P are stable. The lower dashed line is the Maxwell construction with an excess water phase (emulsification failure). G, 
D and P lie on a triple line (compare inset) since they are related by a Bonnet transformation. The right dashed line is the 
Maxwell construction between G and La- Taken from Ref. f|l|. 
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FIG. 8. Phase diagrams as a function of water volume fraction pw = 1 — u and spontaneous curvature co for (a) r = 0.1 and 
n/kBT = 10, (b) r — 0.5 and n/kBT — 10 and (c) r = 0.5 and n/kBT — 2. Two-phase coexistences are indicated by hatched 
regions, (a) has been presented before in Ref. |l|. 
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FIG. 9. Experimental phase diagram for 2:1 lauric acid/dilauroyl phosphatidylcholine and water as a function of water 
volume fraction pw and temperature T (redrawn schematically from Ref. ^). The temperature range from T = 35°C to 
T = 50°C corresponds roughly to the range 0.18 to 0.26 in spontaneous curvature cq in the predicted phase diagram, compare 
Fig. |a. 
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0.794735 
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0.74978 


0.74238 
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0.655994 
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TABLE I. Euler characteristic Xi scaled surface area A* (both in the conventional unit cell) and topology index 
r = {A*'^ /2-k\x\Y^'^ for all TPMS considered, as obtained from local minima of the free-energy functional of a simple 
Ginzburg-Landau theory. Exact values are known for G, D, I-WP, P and C(P) and given in Tab. The structures are 
ordered according to decreasing F; the value for S is probably too large, but difficult to improve numerically. 
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P 


-4 


3A:2 = 2.345103 


0.716346 
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i/k-i = 3.510478 


0.655993 







TABLE II. Scaled surface area A* , Euler characteristic x in the conventional unit cell, and topology index F = {A*'^ /2tt\x\)^^^ 
for those TPMS, for which exact results are available. For surfaces with known Weierstrass representations, we also give A*', 
the number of fundamental domains needed to build up the surface in the conventional unit cell, and a, the resulting lattice 
constant. Here fci — F(\/3/2, \/8/3) where F is the incomplete elliptic integral of the first kind; k2 = K{l/2) / K(\/S/2) where 
K{k) = F{1, k) is the complete elliptic integral of the first kind; fca = K(l/\/2i)/ K{^J2/^)\ and kt, — r(l/3). The surface areas 
and Euler characteristics of P, D and G are related to each other due to their Bonnet-transformations, as described in the text. 
The Bonnet angle 9 = 38.015° for G follows from tan 6 = k2. Note that often a smaller unit cell is chosen for D; then one has 
X = -2 and A* = 1.9188925. 



22 





D 


I-WP 


S 


F-RD 


C(P) 


Sn 


3 83763 


3 46250 


5.41457 
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TABLE III. First moments 5'jv and variance A = S'oS'2/5'1 — 1 of the distributions of Gaussian curvature. In particular, 
So = A* and Si = 2-kx with A* and % given in Tab. If the moments for D are scaled appropriately for a change of unit cell, 
the values given in Ref. are recovered. The moments for G and P are not shown since they follow from those for D due to 
the Bonnet-transformation. The values for A are the same for G, D and P. 
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